home *** CD-ROM | disk | FTP | other *** search
/ Resource for Source: C/C++ / Resource for Source - C-C++.iso / codelib8 / v_10_12 / 1012086a < prev    next >
Encoding:
Text File  |  1995-11-01  |  1.1 KB  |  43 lines
  1. /* _Exp function, Plauger exp() core
  2. corrected treatment of negative
  3.  * arguments -TCP 16.IX.91
  4.  */
  5. #include <math.h>
  6. #include "xmath.h"
  7. #include <float.h>
  8.  
  9. #define hugexp (double)HUGE_EXP
  10. #define invln2 1.4426950408889634074
  11.  
  12. int _Exp(px, eoff)
  13. double *px;
  14. int eoff;
  15. {
  16.     /* Compute e^(*px)*2^eoff, x finite; eoff would be -1
  17.      * for cosh(), sinh() */
  18.  
  19.     double y, g;
  20.     int xexp;
  21.     const static float round[] = {.5, -.5};
  22. #if 0 /* works OK but slow on many machines */
  23.     g = xexp = *px * invln2 + round[*px < 0];
  24. #else
  25.     /* VAX requires -(((short *)px)[_D0]>>15) */
  26.     g = xexp = *px * invln2 + round[((unsigned int *) px)[
  27.         _D0] >> 31];
  28. #endif
  29.         g = (*px - g * (22713. / 32768.)) - g *
  30.             1.428606820309417232e-6;
  31.     y = g * g;
  32. #define c1 15132.70094679458800
  33.     g += g * y * (420.30235984886453 / c1 + y * (1 / c1));
  34. #if FLT_RADIX > 10
  35. #error "rescale for better accuracy"
  36. #endif
  37.     *px = .5 + g / (30265.40189358917686 / c1 - g + y *
  38.         (3362.72154416335329 / c1 + y *
  39.             (30.01511290682109 / c1)));
  40.     return _Dscale(px, xexp + eoff + 1);
  41. }
  42. WRAP_EOF
  43.